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We present a new approach to compute low lying eigenvalues and corresponding eigenvectors for strongly correlated 
many-body systems. The method was inspired by the so-called Automated Multilevel Sub-structuring Method (AMLS). 
Originally, it relies on subdividing the physical space into several regions. In these sub-systems the eigenproblem is 
solved, and the regions are combined in an adequate way. We developed a method to partition the state space of a 
many-particle system in order to apply similar operations on the partitions. The tensorial structure of the Hamiltonian 
of many-body systems make them even more suitable for this approach. 

The method allows to break down the complexity of large many-body systems to the complexity of two spatial sub- 
systems having half the geometric size. Considering the exponential size of the Hilbert space with respect to the geometric 
size this represents a huge advantage. In this work, we present some benchmark computations for the method applied 
to the one-band Hubbard model. 

Keywords: strongly-correlated systems, many-body physics, algorithm, eigensolver, Hubbard model 
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1. Introduction 

In modern solid-state physics strongly correlated many- 
body systems play an increasingly important role, as e.g. 
in the case of the high temperature superconductors, the 
manganites and vanadates, and more recently also in light- 
matter- and ion trap quantum simulators. The numerical 
problem of addressing strong correlations leads to eigen- 
value problems of matrices whose size depend exponen- 
tially on the geometric size and the particle number. 

Several methods have been introduced to solve prob- 
lems of this kind. Quantum Monte-Carlo (QMC, |T,) is a 
powerful method to deal with finite-temperature systems, 
where large system sizes can be reached. On the down- 
side, many system configurations cannot be addressed effi- 
ciently due to the so-called sign-problem [2 ]. Another pow- 
erful method is the Density Matrix Renormalization Group 
(DMRG, [3 ]) which allows to solve for ground-state proper- 
ties of comparatively large-scale electronic structures. Its 
disadvantages include the lack of possibility of treating 
problems apart from ID or pseudo-ID problems. Newer 
methods like the Variational Cluster Perturbation Theory 
(VCPT) also deliver promising results. 

For this work we followed a new approach to break down 
the numerical complexity of such systems. Rather than 
starting from a physical point of view we let ourselves be 
inspired by other fields of numerical simulations, which 



also deal with large-scale eigenvalue calculations. The Au- 
tomated Multilevel Sub- structuring Method (AMLS, [HH]) 
e.g. for linear elastodynamics relies on partitioning the 
physical space into smaller pieces wherein the eigenprob- 
lem is solved as a starting point for the calculation of the 
full system eigenvalues. We refer to our approach as Two 
Sub-system Ground-state Approximation (TSGSA). 

2. Partitioning the occupation number state space 

In second quantization the ab initio and approximation 
free form of the Hamiltonian for the electronic degrees of 
freedom reads 



E 
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where the operator aj, (a„) creates (annihilates) a fermion 
in orbital <f> v . The quantity v = (£, a, a) represents a com- 
bined index including the site- (unit cell-) index i, the in- 
dex a of the basis function to describe the orbitals of the 
various atoms within a unit cell, and the spin a. t and V 
are the corresponding matrix elements for the single par- 
ticle part (hopping) and the interaction, respectively. A 
suitable basis for the many-body problem in second quan- 
tization is the occupation basis where the occupa- 
tion n v for fermions can only be or 1, in accordance with 
the Pauli principle. There are three sources for the intri- 
cacy of a many-body problem. The first one is the number 
of orbital degrees of freedom per unit cell, required for a 
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quantitatively accurate description of the band structure. 
The second concerns the overwhelming number of basis 
states for a true many-body calculation. Let L r be the 
number of unit cells and Lb the number of orbitals within 
a unit cell then there are L = L r Lb different index tuples 
v = (i, a) for each electron spin direction. There are cor- 
respondingly L occupation numbers n v which can either 
be or 1. The number of many-body basis states to de- 
scribe a system of Nf (iVj,) electrons with spin up (down) 
is given by the number of possibilities to distribute JV-j- 
entries 1 and N — TVj- entries among the L occupation 
numbers and likewise for the spin down electrons. Hence 
the number M of many-body basis states is 



M 



L 



L 



It is needless to emphasize that a true many-body ab-initio 
calculation is out of reach. 

In order to study the generic properties of strongly corre- 
lated many-body systems qualitatively it is, however, suf- 
ficient to reduce the number of orbital degrees of freedom 
to a minimum. Common many-body models include up to 
Lb = 3 orbitals per unit cell. 

But even with a strongly reduced number of orbital de- 
grees of freedom there remains a third problem, the struc- 
ture of the interaction part, which is still too complicated 
for an exact treatment, as well by numerical as by ana- 
lytical means. There is reason to believe that the genuine 
many-body effects can already by described and under- 
stood when only short ranged density-density terms of the 
form 



V vv in v h v > 



are retained in the model. The density operator is given 
by n v = a\,a . 

Though not really necessary, but in order to keep the 
number of parameters small, the hopping part is commonly 
approximated by a tight-binding form, allowing for nearest 
neighbor hopping only. The following extended Hubbard- 
model with one orbital degree of freedom per unit cell in- 
cludes an on-site and a next-nearest neighbor interaction: 



Here, we will consider the two most simple and com- 
mon fcrmionic models, the case of spin-less fermions (only 
one spin-species) with nearest neighbor interaction and the 
Hubbard model (only on-site Coulomb interaction). 

For spin-less fermions the Hamiltonian reads 

H = — aja.,- + V ^% ' 
<*j> <*j> 

where i stands for the unit cells and (ij) indicates that the 
unit cells i and j are nearest neighbors. 



The Hubbard Hamiltonian is given by 

Despite of the strong reduction in degrees of freedom, 
the number of many-body basis states is still very large and 
increases exponentially with increasing number of orbitals 
(sites) L. In the case of spin-less fermions we have M = 
(^) and for the Hubbard model, it reads as mentioned 
before M =(£)(£). 

For example, a system of spin-less fermions on L = 20 
sites with N = 10 electrons has M — 184, 756 many-body 
basis states, while a system of half the size, i.e. L = 10 
sites with N = 5 has merely M — 252. For the Hubbard 
model the situation is even more pronounced. Again for 
a half filled system with 20 sites and 10 electrons of each 
species, we find M w 3 • 10 10 , while the half filled system 
on 10 sites has M w 6 • 10 4 . So the numerical complexity 
increases by roughly 10 6 if we double the geometrical size 
of the system. 

The idea of the present paper is to exploit systemati- 
cally the fact that smaller sub-systems have a significantly 
reduced Hilbert space. 

2.1. Sectors in occupation number state space 

The full Hamiltonian, be it the Hubbard or the spin- 
less fermion model, conserves the number of particles per 
spin direction N„. Therefore, the Hamilton matrix is block 
diagonal in the occupation number basis due to the con- 
servation of particle numbers and we solve the eigenvalue 
problem for fixed number of particles (Nf, N±) separately. 
The case of the spin-less fermion model can be traced back 
to the Hubbard case by assuming only one spin direction. 
Here we restrict the discussion to an even number of lat- 
tice sites L = 2Lh and an even number of electrons per 
spin direction N„ = 2^ ;CT . Although the generalization is 
straightforward, it would hamper the discussion unneces- 
sarily. 

Now we split the lattice into two sub-systems A and 
B of equal size. In the approach to be discussed below, 
we will start from the case of completely decoupled sub- 
systems, i.e. all hopping and interaction effects between 
the two sub-systems are ignored. Consequently, the num- 
bers of electrons per spin direction are conserved in each 
sub-system separately. Due to the particle number con- 
servation in the sub-systems, the corresponding occupa- 
tion number state space can be split into different sectors 
(If, li), which are characterized by the number of particles 
per spin direction (Nf = N h j + Z t , A^j. = N h> ± + Zj.) in 
the sub-systems under consideration. Here we have intro- 
duced a the quantity l a which specifies the deviation of the 
actual particle number for spin a from the reference value 
Nh, a - Obviously, l a can range from —Nh !<7 to +Nf l , (J . 

Next we construct a complete basis for the entire sys- 
tem, by tensor-products of the eigenvectors of the decou- 
pled sub-systems. In order to achieve the correct parti- 
cle numbers per spin direction for the entire system, only 
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Figure 1: Particle number sectors £4.)) for a Hubbard system. 
Shells of sectors with the same Manhattan have the same color. The 
yellow line marks a possible combination of sectors to build a parti- 
tion. 



specific sectors (l a ,f,l a i) of the two sub-systems can be 
combined, (If, Zj,) for sub-system A and (—If, — Zj.) for sub- 
system B. It will turn out expedient to introduce the 
Manhattan distance in the (l a ,f , Z a ^)-plain, which is given 
by d((l a ,t,lai) = \la,t\ + l'a,4.|- Sectors of the same Man- 
hattan distance form shells, which are of equal importance 
as far as eigenvectors and eigenvalues are concerned. 

In fig. [T] different shells up to to Manhattan distances 
2 are depicted. As pointed out before, within one shell 
there are always two opposite (If, Z^)-pairs, one for each 
sub-system, which form a particle number partition which 
is used in the tensor-product basis for the full system. One 
such example is indicated in the figure by the yellow line. 
The central square denotes the sector with Manhattan dis- 



tance 0, i.e. (If = li = 0). 

2.2. Hamilton operator of the sub-systems 

It is expedient to introduce a combined site index i 
(a, i), where a G {A, B} refers to the sub-systems and and 
henceforth i G {l,...,Nh} enumerates the sites within 
each sub-system. The Hamiltonian can be decomposed 
into three parts 

H = H A + H B +H AB . 

■=H 

For spin-less fermions the parts of the Hamiltonian are 



a e {A, B} 



Hab = 



:H" 



where indicates that the indices i,j have to be chosen 
such that (A, i) and (B, j) belong to nearest neighbor sites. 
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Figure 2: Partitioning the occupation number state space by num- 
ber of particles in each sub-system. Here, representatives of each 
partition are depicted. 



For the Hubbard model the splitting yields 
H a = -t +U^2n atit h aM ; a G {A,B} 



Hab = -t^2 i^A.i^Bj, 



B,j,u u A,i,(T 



nA,inB,3 



In case of the spin-less fermion model there is a two- 
particle coupling term between the two sub-systems (H™ B ) 
which is detrimental for cluster perturbation theory (CPT 
or VCA, [3 [7])- In the present approach it does not make 
any difference at all. 

Now, the particle numbers in the two sub-systems are 
not conserved, due to the inter-sub-system hopping, and 
in principle all particle numbers between and TV are con- 
ceivable for each sub-system. However, the most probable 
number of particles in the sub-systems is Na = Nb = Nh ■ 
The particle number fluctuations in the two sub-systems 
is subject to the constraint Na + Nb = N and likewise 
for the two spin species Na,o- + N B .a = N h ^ a . If the two 
sub-systems are decoupled, i.e. for Hab = the particle 
numbers in each sub-system are conserved as well. 

In the case of the Hubbard model the partitioning con- 
cerns both spin species separately. A partition (Zfjij.) 
therefore represents the situation of Na.o — Nha + la and 
Nb,<j = Nfio ~ ten respectively. Typical configurations for 
partitions with I = and I = ±1 are depicted in fig. [2j 

The eigenvectors of H are simply tensor products of the 
eigenvectors of the separate sub-system Hamiltonians H a 
and the eigenvalues are given by the sum of the correspond- 
ing eigenvalues. Let the eigenvalue problem of sub-system 
a with tV q = N + I electrons be given by 

H«\*P) a =EP\*P) a . 

Here the meaning of the the indices is as follows. The lower 
outer index a stands for the sub-system, the upper index 
(I) represents the particle number (N a = Nh + 1), and the 
lower index v enumerates the eigenvalues and eigenvectors 



3 



for this particle number. Note, that the eigenvalue spec- 
trum is the same for both sub-systems if they are occupied 
by the same number of electrons. 

The eigenvectors and eigenvalues of Hq in partition I are 
given by 



2#> := E P A + E(-J), 



where the index v represents the entire system with v = 
[y Al v B ). 

Due to the missing coupling of the sub-systems, the 
eigenvalue problem of Hq has a strongly reduced complex- 
ity. If we add Hab by perturbation theory, there will be 
a first order contribution only from the interaction term 
in the spin-less fermion case. The hopping terms do not 
contribute in first order, as they change the number of 
particles in the sub-systems. 

The ground-state of H is obtained for Na — Nb = N^, 
i.e. it belongs to partition I = 0. The second order energy 
correction is given by 



N h 
l = -N h 



\(4 1) \Hab\^ 0) ) 



E 



(1=0) 



E, 



As usual indicates that (I = OAv = 0) is excluded from 
the sum. Only I ~ and / = ±1 contribute to the energy 
correction, because only one electron can hop at a time 
across the border . In other words, the partitions I = ±1 
come into play. The importance of the unperturbed eigen- 
vectors of partitions I = ±1 is determined by the matrix 
the element M^^ :— (^^^^Hab^q^) and the in- 
verse of the energetic distance E^~ ^ — E^ . The latter 
implies that low eigenstates of the unperturbed system are 
more important for the ground-state of the entire system 
than those with higher energies. In addition the kinetic 
coupling matrix element is driven by the hopping across 
the border and yields for the spin-less fermion model 

ij 
ij 

Hence the relevance of eigenvectors of the sub-systems 
for the partition I = +1 depends on the occupation of the 
border sites in those states. We obtain a similar result 
for Mu~ . From these considerations we can derive a 
criterion for the importance of the contribution of excited 
eigenvectors of the sub-systems to the total ground-state 
of the entire system. 

We observe that the partitions I — ±1 contribute to 
second order energy or rather first order vector correction. 
The second order correction for the eigenvectors depend on 
the partitions |^| = 2 and generally the n th order correction 
requires |^| = n. This conclusion can also be obtained by 
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Figure 3: Block-structure of the Hamiltonian for the spin-less fermion 
model. As explained in the text, Na = + I and Nb = N^ — I. 
Only partitions up to I = ±2 are depicted. White spaces represent 
blocks with zero matrix elements. 



starting a Lanczos iteration with the eigenvectors of the 
unperturbed system. Each Lanczos iteration then increase 
the required partition |/| by one. 

In the Hubbard model we have to deal with the pair 
I = Starting from I — (0,0), to which the ground- 

state of Ho belongs, each Lanczos iteration modifies one 
of the values l a by ±1. I.e. the first iteration includes 
all partitions with Manhattan distance 1 from the center 
/ = (0,0), the second iteration adds all partitions with 
Manhattan distance 2 and so forth. 

3. Block structure of the Hamiltonian 

The partitioning described in the previous section leads 
to a natural block structure of the Hamiltonian, which is 
depicted for the spin-less fermion model in fig. [3] 

The upper-left block corresponds to partition I = that 
contains all many-body basis vectors in occupation num- 
ber representation \{n}} with equal number of electrons in 
both sub-systems. 

In general the basis vectors in partition I are a tensor 
product of the contributions of the sub-systems 



where v — (va,vb) enumerates the basis vectors of the 
entire system, while v a enumerates the basis vectors of 
sub-system a. The occupation number basis vectors of the 
two sub-systems have an outer index a £ {A, B} indicating 
the sub-system they belong to and they are constructed by 
the corresponding creation operators 



i=l 



y,n v (a,i) =Nh±l; n v (a, i) £ {0, 1} , 

i 

where |0) represents the vacuum vector of sub-system 
a. As a consequence of the representation by creation 
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Figure 4: Block-structure of the Hubbard Hamiltonian. The mean- 
ing of l a is explained in the text. Partitions up to Manhattan dis- 
tance 1 are depicted. White spaces represent blocks with zero matrix 
elements. 

operators, the two factors of the tensor product do not 
commute. Commuting the factors yield an additional sign 
(— X) A B . Similarly, there could be an additional sign 
when computing matrix elements in the tensor basis. Due 
to the tensor structure of the basis the computation of 
the matrix elements of H can be simplified significantly. 
The contribution of H to the diagonal block of partition 
Z reads 

{Ho)^ = {n^\H A \n^) A 5 v . B , UB 

+ ^A^ ) \H B \n^) B . 

As the diagonal blocks contain a fixed number of elec- 
trons in each sub-system there is no contribution stemming 
from the hopping part of Hab- However, in the spin-less 
fcrmion case the interaction term H^b contributes to the 
diagonal block as well 

(hTbYH = V Tl «\^,i\n { Jl) A (nl-»\n B M- B l) ) B • 

ij 

In both models under consideration the off-diagonal 
blocks are solely due to the hopping part in Hab, which 
changes the number of particles in the sub-systems by ±1. 
Consequently, only those blocks possess non-zero entries 
for which the partition indices differ by \V — l\ = 1 in the 
spin-less fcrmion case. In the case of the Hubbard model 
the condition for non-zero blocks reads |Z^ — l^\ + \l± — 1'^\ = 
1. I.e. in the Z-j-, Zj, plain only those sites are coupled which 
have Manhattan distance 1. The corresponding spin re- 
solved block structure is depicted in fig. [4j 

For the Hubbard model the matrix elements in the off- 
diagonal blocks with Z' > Z are 

Due to the reordering of the creation/annihilation oper- 
ators there is an additional sign factor s = (— 1) h+ . 



The matrix elements for Z' < Z follow from the hermiticity 
= (M^ 1 , jP) . The matrix elements for the spin-less 
fermion case are obtained by restricting the sum over the 
spins to one spin species, to a =t say. 
All block matrices have the structure 

mW=J2(m^ a (m^ b , (1) 

K 

which will be exploited later on. 

4. Ground-state eigenvalue problem 

In the present approach we first solve the eigenvalue 
problem of the decoupled sub-systems for the different par- 
titions. 

(H ) (1 ' 1) VW =V^DW , 

where is the unitary matrix of eigenvectors and the 
diagonal matrix D™' contains the eigenvalues for partition 
Z. As before, in the case of the Hubbard model Z stands 
for Z = (Z t ,Zj.). 

We begin with Z = (Z-j- = Z^ = 0) and include gradually 
partitions of increasing Manhattan distance. We always 
include all partitions to a a given Manhattan distance, 
like (l f , Z4.) € {(1, 0), (-1, 0), (0, 1), (0,-1)}. From the set 
of eigenvectors { | V^)} of these partitions we keep a certain 
number, which we call cropping number. Next we form the 
eigenvector of the entire system by a linear combination 
of those selected eigenvectors of the partitions included so 
far, i.e. 

i^hEE^i^) ■ (2) 

1=0 v 

If all partitions and all corresponding eigenvectors are in- 
cluded the result will be exact. Below, we will give a de- 
tailed study of the convergence as far as the upper parti- 
tion number Z* and the cropping numbers are concerned. 

Given Z* and the set of retained eigenvectors {{tpi^}} 
for Z < Z*, the coefficient vectors cfj is given by the 
eigenvectors of the matrix 

H ( Si = {4' ] \H\^) , (3) 

spanned by the restricted set of unperturbed eigenvectors 
{IV^)} f° r I < I* ■ The approximated matrix H has the 
same block structure in the partition numbers as the orig- 
inal matrix. Moreover, the matrix ,l > of block Z', Z is 
given by 

H (i',i) = (y(i')^ H {i' -i)y<,i) ) 

where the matrix V"^ contains column-wise those the or- 
thonormal eigenvectors of the diagonal block -Hi ? , which 
are used in the expansion (12]). 




Figure 5: Illustration of the construction of the TSGSA basis V. 
Each partition corresponds to a fixed number of spin-up and spin- 
down particles in both sub-systems, which eigenvectors are combined 
tensorially (blue) to Vw. For each sub-system A and B the vectors 
are obtained by solving the corresponding eigenvalue problem. The 
eigenvectors are truncated e.g. by keeping only the vectors on the 
lower end of the spectrum (red) . 



The contribution of Hq to the diagonal block ffO'O is 
diagonal, containing the retained eigenvalues D"'. For the 
computation of the block matrices resulting from Hab one 
can exploit the tensor structure outlined in eq.Q along 
with the tensor structure of 

y{l) = y(A,l) y(B,-l) 
V U',V V u'., VA V' n ,VB ' 
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0.0019 



Table 1: Ground-state energies of a Hubbard system (L = 12, jV-j- 



N, 



bc.,U = 10) as estimated by TSGSA compared to Exact 



Diagonalization (Lanczos). The exact solutions is Eo = —3.232383 
and the basis size of the full Hamiltonian is 853,776. The first column 
shows the croppings for the individual partitions used (e.g. for the 
first line 80 eigenmodes were used for the primary partition 1 = 0; 
60 for the partition I = ±1, and so on). 
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N a 


Croppings 


TSGSA 


DMRG 


rel. error 


4 


2 


4/2 


-0.882 


-0.911 


0.032 


8 


4 


36/24 


-1.937 


-1.975 


0.019 


12 


6 


50/25 


-2.949 


-3.041 


0.030 


16 


8 


50/25 


-3.961 


-4.109 


0.036 


20 


10 


50/25 


-5.021 


-5.178 


0.030 


24 


12 


50/25 


-6.090 


-6.245 


0.025 



Table 2: TSGSA for multiple Hubbard system sizes (obc), compared 
with DMRG. 



eigenvectors 



V = 






V i 
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where the columns of V a ' 1 contain the eigenvectors of H a 
retained by the cropping process. So we see that all op- 
erations are restricted to vectors and matrices of the size 
given by the sub-systems. 

Fig. [5] illustrates the construction of the overall unitary 
matrix V. 

The overall unitary transformation from the orthonor- 
mal occupation number basis to the orthonormal basis of 
the eigenvectors of Hq has a block structure, correspond- 
ing to the partitions (/<[-, Zj,). The partitions are enumer- 
ated with increasing Manhattan distance and within a 
shell of fixed distance d counter clockwise beginning with 
(l f = d,k = 0). 



V 






V : 











^(3) 



■7 



where the columns are the eigenvectors of Hq in the vari- 
ous partitions. Zeros are the corresponding zero matrices. 
If we only retain the restricted set of eigenvalues we ob- 
tain a Nf, x m matrix, where m is the number of retained 



the zero matrices are adjusted appropriately. V is no 
longer unitary, but V^V — I still holds for / being the 
m x m identity matrix. The matrix P := VV* is the 
projection matrix into the space spanned by the retained 
eigenvectors. The above procedure corresponds to the 
eigenvalue problem of the projected Hamiltonian matrix 
H = PHP. 

5. Numerical results 

In order to determine how well the ground-state of a 
Hubbard model is approximated by the present approach 
a series of numerical simulations was performed and com- 
pared to exact values achieved by standard algorithms like 
the Lanczos method. 

As can be seen in tab. [T] the exact eigenvalues are 
approximated well by a comparably small vector space. 
Furthermore, the approximations can be enhanced by us- 
ing a limited amount of additional eigenmodes of the sub- 
systems. 

A study for different system sizes can be seen in tab. 
[2] Here, the DMRG implementation by Reinhard Noak 
was used for comparison. Note, that the relative error is 
stable to decreasing with system size for constant number 
of eigenmodes taken into account. This is due to the fact 
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t[s) 


TSGSA 


rel. error 


25/25/25 


5 


2223 


-5.5101 


0.066 


50/50/50 


10 


5123 


-5.5305 


0.061 


75/75/75 


10 


15916 


-5.5322 


0.061 


100/100/100 


20 


23369 


-5.8169 


0.013 



Table 3: Convergence of the method for a Hubbard-type system (L = 
22, N-^ = N l = 11, pbc, U = 10), compared with DMRG (ground- 
state energy: —5.8907). N COICS indicates the number of computation 
cores used, t gives the total run-time in seconds. 



that with increasing system sizes the primary partition be- 
comes more important and hoppings over the sub-system 
boundary have less weight. 

The largest Hubbard-type eigensystem calculation by 
exact diagonalization known to the author was preformed 
at the Earth Simulator j8j. The Hamiltonian used (L = 22, 
N t = 9, JVj. = 8) had 1.59 • 10 11 unknowns. Tab. [| shows 
the convergence of a system (L = 22, Nf = N± = 11) 
for the for physics important case of half-filling evaluated 
by the presented approach. The corresponding full basis 
size exceeds the mentioned world-record size by a factor 
of 3 (5.0 • 10 11 ). The calculations were performed using 
a parallelized implementation on a 8 QuadCorc-Opteron 
CPU cluster at Graz University of Technology. 

6. Implementation 

A key to the success of the algorithm is to find a scheme 
for calculation of more than just a few lowest eigenvec- 
tors, say in the order of 100. For this purpose an Im- 
plicitly Restarted Lanczos algorithm was used. Note, that 
the final, transformed Hamiltonian does not have to be 
calculated explicitly but a method for its application on 
a vector can be derived using the sub-system eigenvalues 
and eigenvectors. 

Before combining the individual partitions to the effec- 
tive Hamiltonian the calculations of the individual sub- 
systems can be done completely independent. So a paral- 
lclization of the calculation scheme could be done rather 
easily by dividing the different sub-system occupation con- 
figurations to different computation cores and solve the 
according eigensystems without need for communication. 
For larger systems this may not be practical, as the imbal- 
ance of numerical complexities among the different sub- 
systems makes a straight forward parallelization more in- 
efficient. For this cases a parallelized version of the Im- 
plicitly Restarted Lanczos using algorithm was used. 

The algorithm was implemented in the C++ program- 
ming language, as a parallelization framework OpenMPI 
was used. 

7. Conclusions and Outlook 

The presented approach may lead to a new way of per- 
forming calculations in strongly correlated material sci- 
ences. The results are promising compared to exact di- 



agonalization. Although other sophisticated methods like 
DMRG, VCPT, or QMC exist, the TSGSA has the advan- 
tage of producing an explicit matrix representation. Fur- 
thermore it is not limited to one dimensional systems. It 
can be easily extended to more dimensions. In this case 
the off-equilibrium partitions may become more important 
due to larger interfaces between the sub-systems. 

A further possibility for developing the algorithm are an 
intelligent way of selecting eigenmodes of the sub-system. 
It can be shown that many of them do not contribute to 
the full system ground-state. 
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